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Abstract. We have developed Monte Carlo wave function simulation schemes to study cold collisions 
between magnesium atoms in a strong red-detuned laser field. In order to address the strong-field problem, 
we extend the Monte Carlo wave function framework to include the partial wave structure of the three- 
dimensional system. The average heating rate due to radiative collisions is calculated with two different 
simulation schemes which are described in detail. We show that the results of the two methods agree and 
give estimates for the radiative collision heating rate for 24 Mg atoms in a magneto-optical trap based on 
the So- Pi atomic laser cooling transition. 

PACS. 32.80.Lg Mechanical effects of light on atoms, molecules, and ions - 32.80.Pj Optical cooling 
of atoms; trapping - 42.50.Lc Quantum fluctuations, quantum noise, and quantum jumps - 2.70.Uu 
Applications of Monte Carlo methods 
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1 Introduction 

Trapping and cooling of neutral atoms in magneto-optical 
traps has become an important tool in atomic physics for 
spectroscopic and fundamental studies since it was first in- 
troduced in late 1980's pQ. In such traps atoms are cooled 
down to the milliKelvin region and below it (depending 
on the particular element), and subsequent evaporative 
cooling in magnetic traps has made it possible to reach 
the quantum degeneracy limit, where bosonic atoms tend 
to form Bose-Einstein condensates [2] , and fcrmionic ones 
demonstrate the exclusion principle by filling the Fermi 
sea |3I4) . Despite the success in going colder and denser, 
the behavior of the atomic cloud in basic magneto-optical 
traps still remains only a partially understood problem. 
One of the challenges is to understand if and how the 
atom-atom interactions in the form of inelastic cold atomic 
collisions affect the thermodynamics of the trapped and 
cooled cloud of atoms. Studies of such collisions have, 
however, also a more fundamental role, because cold col- 
lisions in the presence of light challenge the traditional 
approaches of collision and scattering theory |5l(:i| . 

In general laser cooled and trapped neutral atoms pro- 
vide an excellent system to study the basic aspects of col- 
lision theory. The slow relative motion of atoms guarantee 
the validity of the Born-Oppcnhcimcr approximation and 
one can consider the collision as a one-cycle vibration of 
a quasimolecule formed by the colliding atoms (we prefer 
the term "quasimolecule" to "molecule" because any bind- 
ing of atoms to each other is only temporary). One can 
separate the relative motion of two atoms into radial and 



angular parts, and the number of involved partial waves 
is limited to only a few lowest ones. Especially for the 
main isotopes of the alkaline earth atoms, such as Mg, Sr 
and Ca, the missing nuclear spin leads to a single, non- 
degenerate atomic ground state, and the bosonic nature 
of these atoms forces an even parity to the partial waves 
(or, in the quasimolecule picture, to the ground state rota- 
tional quantum numbers). Naturally, in the limit of very 
low temperatures only the s-wave is allowed, and under 
suitable conditions the atoms can form a Bose-Einstein 
condensate [7]. 

The simple picture given above becomes distorted if we 
add the presence of light. Firstly, on a single atom level 
the absorption of a photon, although needed for laser cool- 
ing, will eventually limit the temperatures that one can 
achieve (the recoil limit) [H] . Thus ultracold temperatures 
can be obtained either by evaporative cooling or by us- 
ing very sophisticated absorption-free techniques such as 
velocity-selective coherent population trapping [H], or Ra- 
man cooling [TP] , Secondly, the laser can also excite the 
quasimolecule formed by the colliding atoms [11112113] . 
The result of this process depends very much on the fre- 
quency and intensity of the laser light. In the basic laser 
cooling situation, however, the light field is usually near- 
resonant and strong. Unfortunately, this is also the region 
where theoretical treatment of the problem is most diffi- 
cult, and where the assigning of experimentally observable 
signals to specific processes is hard. 

The main bulk effects arising from collisions in the 
standard laser cooling and trapping process in a magneto- 



2 



J. Piilo et al.: Cold collisions in strong laser fields: partial wave analysis of magnesium collisions 



optical trap are assumed to be a) the loss of atoms from 
the trap, and b) the heating of the atomic cloud. The first 
process has been studied in detail both experimentally 
and theoretically, but mostly at weak fields (in a typical 
experiment the cooling lasers are shut off and replaced 
by a highly controllable, weak probe laser); for reviews, 
see Refs. [5I6| . The second process has been almost com- 
pletely ignored due to the difficulty of modelling it |14) . 
and to the difficulty of separating the collisional heating 
from other heating effects such as reabsorption of scat- 
tered photons jT5| . 

In this paper we build a dynamical description of colli- 
sional heating by a strong laser field in a three-dimensional 
two-atom system. Our tool is wave packet dynamics |lfi| . 
i.e., time-dependent description of the collision process. 
Spontaneous emission of photons in the semiclassical laser 
field approximation makes the two-atom system an open 
quantum system |17j . and we use the Monte Carlo quan- 
tum jump method as a tool ^B]- We take the dominant 
magnesium isotope 24 Mg as our study case, because the 
lack of hyperfine structure reduces the number of molecu- 
lar states and thus allows a quantitative study ^Hj , where 
the complexity arises only from the partial waves. Also, 
laser cooling and trapping of magnesium has been achieved 
experimentally |2()I21| . 

Previously near-resonant collision studies have been 
limited to two-state models, and at best the partial waves 
have been treated independently (ignoring second-order 
couplings) but at strong fields this is not possible as 
several waves are coupled to each other. We have recently 
taken a step to the direction of including the partial wave 
structure in a more realistic way to the description of cold 
collisions in strong fields [22] - Of course, it is not possi- 
ble to include all partial waves in a realistic study, so the 
partial wave manifold must be truncated. 

In building the framework for our wave packet calcu- 
lations we have noticed that the artificial truncation of 
the infinite sequence of partial waves must be done with 
care, otherwise it may lead to significant errors. Another 
problem in cold collisions is the setting of the initial state 
for the atom-atom system, as the quasimolecule states are 
strongly coupled even asymptotically. After outlining the 
background for inelastic processes and the properties of 
the magnesium system in Sec.|21 we discuss these problems 
and our solutions to them. In Sec.|3we discuss how to con- 
nect our simulations with a heating rate. By performing 
the wave packet calculations in our framework, we show 
that it is possible to assign a heating rate to the system. 
We also study the intensity-dependence of this rate. These 
results are reported in Sec. 01 where we also show that the 
long-time behavior of the atom-atom system, i.e., accumu- 
lation of relative kinetic energy in several subsequent col- 
lisions, can be understood and predicted reasonably well 
from the short-time one-collision data. Finally, in Sec. 
we discuss our results and their consequences. 
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Fig. 1. The quasimolecule view of the colliding magnesium 
atoms in a laser field. 

2 Framework 

2.1 Inelastic processes 

Let us consider the situation in Fig. Q where we show 
schematically the basic molecular potentials of our model, 
as a function of the atomic distance R. For best effi- 
ciency the frequency of a cooling laser (u>l) is usually set a 
few atomic electronic transition linewidths (-T a t in energy 
units) below the exact resonance (u)q), i.e., detuned to the 
red side of the transition [detuning S = h(u>L — loq) < 0]. 
A single atom is thus slightly off-resonant with the laser 
light. However, we can also view the situation in the quasi- 
molecule picture. Then one of the molecular states can at 
some specific atomic distance Rc (Condon point) couple 
resonantly to the molecular ground state. Usually there 
are many such states, but for 24 Mg (and also for the other 
major isotopes of alkaline earth atoms) there are only two, 
namely 1 S+ and 1 II g |1HI23I24I25] . Since the excitation of 
the 1 n g state is strongly forbidden at long distances, the 
system reduces in practice to a simple system consisting of 
two electronic states. Note that asymptotically the molec- 
ular ground state ( 1 S+ here) corresponds to both atoms 
being in the atomic ground state, and the excited molec- 
ular state corresponds to one atom being in the atomic 
ground state (here 1 So) and the other one in the atomic 
excited state (here 1 Pi). 

Eventually at very short atomic distances the poten- 
tial for the excited molecular state becomes repulsive, and 
thus below the dissociation limit we have a set of bound 
vibrational states (dotted lines in Fig. ^| . When the cool- 
ing laser is just below the excited state dissociation limit, 
the vibrational line broadening is usually too large for the 
states to be resolved. By increasing the laser detuning, 
the states become resolved and one can perform spectro- 
scopic studies of these states. This leads to photoassocia- 
tion spectroscopy [Fd], which is a very useful tool for deter- 
mining the s-wave scattering lengths for atomic collisions, 
and for creating translationally cold molecules from cold 
atoms. 

If we describe the situation in the dynamical picture 
based on only the electronic molecular states, the inability 
to resolve the vibrational states means simply that if the 
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quasimolecule is excited at Rc, it has a high probability to 
decay back to the ground state via spontaneous emission 
of a photon, before reaching the small values of R. If the 
quasimolecule survives to the region of small R, its state 
may change into a one with lower energy (asymptotically) , 
and the energy difference between these two excited states 
goes into an increase of the relative kinetic energy. This 
increase is large enough for the atoms to be lost from 
the trap, and the mechanism is called the fine-structure 
change (FS) for alkali metal atoms jl 1112118] . and state 
change (SC) for alkaline earth atoms |19| . 

Another possibility is that the quasimolecule decays 
before reaching the small R region, but has acquired enough 
relative kinetic energy for the atoms to exceed the trap 
depth (usually on the order of 1 K). In magneto-optical 
traps Doppler cooling works efficiently on hot atoms, but 
if the trap depth is exceeded, even this cooling can not 
recapture the hot atoms produced in collisions. Both ef- 
fects lead to a non-exponential decay in the number of 
trapped atoms and can be observed via monitoring the 
number of atoms (collisions with background atoms etc. 
lead to exponential decay). This is called radiative escape 

(RE) EH]. 

Radiative heating, which is studied in this paper, fol- 
lows the radiative escape mechanism, but occurs for atoms 
that have not managed to exceed the trap depth in their 
kinetic energy increase before decay. The loss mechanisms 
are very visible and important, but compared to them the 
radiative heating is really a bulk effect, especially at strong 
fields. The resonant excitation at Rc is very strong, but 
the decay also concentrates to its vicinity. At strong fields 
the population that decays near Rc can also undergo re- 
excitation |14|. i.e., we can not separate the dynamics into 
a sequence of separated excitation and decay processes, as 
has been done in the weak field studies. The inclusion of 
partial waves complicates the situation further, as we shall 
discuss in Sec. 12.51 
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Fig. 2. Partial wave I and rotational state J coupling scheme 
for the ground and excited states. 



where k r is the wavenumber of the 1 So- 1 Pi transition, R 
is the internuclear distance, and -T at is the corresponding 
atomic linewidth (r at /h — 78.8 MHz). For short range we 
actually flatten the excited state potential when it begins 
to correspond to kinetic increases that exceed the max- 
imum values allowed by numerics (and survival to short 
range is small in any case) . The position-dependent molec- 
ular linewidth including the relativistic retardation correc- 
tions is 

r{R) = T at jl - — ^ 3 [(k r R) cos (k r R) - sin {k r R)\ j . 

The wave number k r forms the basis of the recoil unit 
system, which is used in the numerical computations. The 
recoil energy is E r = h 2 k 2 /2fi, where fi — iriM g /2 is the 
reduced mass of the two-atom system, and the recoil unit 
of length is \/2ir — k" 1 , where A is the transition wave- 
length corresponding to Tiloq. For the 1 So- 1 Pi transition 
of 24 Mg, E R = k B x 9.8 fiK and A = 285.21 nm Q]5]. 

Within the rotating wave approximation we can shift 
the ground state up in energy by the photon energy Hlul, 
giving the atomic states (and asymptotically for the molec- 
ular states) an energy difference equal to the detuning 
5. Thus the Hamiltonian of the system, with the matrix 
structure due to the partial waves made explicit, reads 



2.2 Magnesium quasimolecule 

We have selected as our study case the 1 So- 1 Pi laser cool- 
ing transition of 24 Mg. The methods and qualitative con- 
clusions apply to other alkaline earth atoms, and to a 
lesser degree also to other elements where, however, hyper- 
fine structure will complicate the picture. For the system 
of two 24 Mg atoms, the ground electronic quasimolecule 
state is X S+ and, as discussed above, for excited state we 
choose the state, for which the P and R rotational 
state transition branches are possible, see Fig. [2 Since for 
symmetry reasons only even partial waves exist for the 
ground state, only the odd ones contribute in the excited 
state due to the missing Q branch. 

For small detunings the Condon point Rc lies in the 
long range where the ground state is approximately flat, 
and the excited potential is ^5] 



U(R) 



3r at 

2 {k r Rf 



[cos (k r R) + (k r R) sin (k r R)] , (1) 



H 



, + Ue.l-l 

Vi,i-x 




t/ki 



Vi-i.i 

+ Ug.l 

Vi+i,i 





Vi,i+i 
n + £4,z+i 



where [/ki n is the kinetic energy 



h 2 d 2 

2fi dR 2 ' 



(3) 



(4) 



and the potentials U g and U e include the atomic interac- 
tion potential U{R) from Eq. Q and the centrifugal part: 



U, 



9,1 



h 2 i(i + i) 

2nR 2 : 



U(R) 



n 2 j(j + 1) 

2^i? 2 



(5) 
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Table 1. Factors a.ji for the couplings between the various 
partial waves/rotational states. I denotes the ground state par- 
tial wave angular momentum number and j the corresponding 
excited rotational state. P denotes the branch where the par- 
tial wave I couples to the excited rotational state with angular 
momentum j = I — 1, and R denotes the branch coupling I to 
j = l + l. 



branch 


ajl (I = 0) 


a* (I + 0) 


P 
R 


v / 2/ ? 3 


/ 1/3 
V 2l + 1 
/(i+l)/3 
V 2i +! 



The couplings between the ground state partial waves (I) 
and excited state rotational states (j) are denned as 

^=^y^£, (6) 

where the factors ctji are given in Table ^ and fl is the 
Rabi coupling in energy units. Following Ref. |19| . we can 
relate the Rabi coupling to laser intensity as 

= 0.5304 ^I{W/cm 2 ), (7) 

^ at 

where / is the laser intensity. The previous weak field 
work m3 was performed at the intensity 1 mW/cm 2 , where- 
as our work concentrates on values £2/r at ~ 1, indicating 
the intensity region 1 W/cm 2 , i.e., about three orders of 
magnitude stronger fields. 

2.3 Monte Carlo simulations 

The Hamiltonian © is not sufficient for describing colli- 
sions in the presence of laser light, because it does not take 
into account the spontaneous decay of the excited state. 
At weak fields it is unlikely that the decayed population is 
re-excited, so one can add decay to the Hamiltonian as 
an imaginary potential iT/2 for the excited quasimolecule 
states, with the position-dependent F given by Eq. (0). 
The resulting problem can be solved with traditional time- 
independent scattering methods for trap loss |14I19| . For 
radiative heating we need to estimate also the kinetic en- 
ergy change that took place before the decay, and in strong 
fields we also need to allow for re-excitation. This forces us 
to treat the problem as a time-dependent one, with initial 
and final states, rather than with boundary conditions for 
ingoing and outcoming single-energy quantum waves. In 
this picture, the colliding atoms are described by a wave 
packet of relative motion, initially peaked around some 
mean value of relative velocity v. When we refer to the 
initial velocity of the wave packet, it actually means (v); 
the same applies to the kinetic energy E^. 

Since the wave packet system is now an open quan- 
tum system, its description would require one to use a 
density matrix instead of a wave function. That would 
also double the spatial degrees of freedom [density matrix 



p(R,R';t) vs. wave function <P(R,t)), taking the problem 
beyond current computational capacities [lfi| . However, 
we can nevertheless use the wave function description. In 
the Monte Carlo wave function scheme ^B], the decay is 
incorporated in a stochastic manner. We add the iT/2 
term to the Hamiltonian In addition, at each time 
step of duration 6t in the numerical propagation of the 
wave packet with time-dependent Schrodinger equation, 
the wave packet is allowed to collapse onto the ground 
state with a probability 

Pjump = stJ2[ dR{Vi(R)}*Vi(R)r(R)/7i, (8) 

i 

where the sum goes over all excited quasimolecule states 
(here all the rotational states j). The occurrence of the col- 
lapse (a quantum jump) is determined by a random num- 
ber. At the collapse \P g (R,t) is deleted, \P e (R,t) becomes 
the new ^ g (R, t) , and then the excited state component is 
set to zero |27| . Whether a jump actually occurred or not, 
the normalization of the wave function is restored at each 
time step. 

The system properties are calculated as ensemble av- 
erages over several such random realizations. The benefit 
of the strong fields regime is that the excited state occu- 
pation is high, and jumps are frequent; this improves the 
statistical accuracy and the relevant properties, including 
the average kinetic energy as a function of time, can be 
obtained with reasonably small ensemble sizes. Typically 
our ensembles have 64 or 128 members, which is quite suf- 
ficient for calculating the bulk properties such as average 
kinetic energy |14I26| . In cold collision studies this method 
has been applied previously to two-state models |I4l2fi| . 
to collisions in optical lattices 28 29 , and an initial step 
towards including the partial waves to cold collision prob- 
lems was presented in Ref. |22| . 



2.4 Semiclassical picture 

The wave packet approach treats the molecular dynam- 
ics quantum mechanically. The laser field is assumed to 
be classical, though, so it appears as a dipole-mediated 
coupling of the electronic states. Decay via spontaneous 
emission is included with the quantum jump approach. 
Thus the reabsorption of a spontaneously emitted photon 
is ignored, but this does not really play a role in colli- 
sion dynamics. Also, at the temperatures considered here, 
we can also ignore the photon recoil for absorption and 
emission events. 

There is, however, an aspect of classical dynamics that 
one can add to the system. Initially the wave packet is 
rather well defined in position and momentum, so it is 
possible to define a classical collision trajectory for the 
approaching atoms. In other words, the excitation can be 
treated as a dynamical process where the quasimolecule 
traverses a level crossing at Rc along a classical trajec- 
tory. If only two states are involved, and decay is ignored 
for the brief moment of passing Rc, one can apply the 
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Landau-Zener model 5 14 . This model predicts an exci- 
tation probability 

^-"(-s^b)- (9) 

where vc is the relative momentum, Vc is the coupling 
between the molecular states, and U'(Rc) is the slope of 
the potential difference, all calculated at the crossing point 
Re 

Even with its limitations, this simple model enables 
one to define the concept of a " strong field" . In the atomic 
case, a field is strong if one is at the saturation limit of 
the atomic transition. The on-resonance saturation inten- 
sity is usually defined as 1 I s — 7rft,c/(3A 3 r), where r 
is h/r at . For the magnesium 1 So- 1 Pi transition one gets 
I s = 0.444 W/cm 2 , which is of the same order as the inten- 
sities used in our simulations. But one should note that the 
quasimolecular excitation probability can be close to unity 
even when the atomic excitation is far from being satu- 
rated. The Landau-Zener expression indicates that slow 
motion or gentle potential slope at Rc can compensate in 
situations where the ratio V/5 is otherwise small. There- 
fore the concept " strong field" must always be conditional 
to temperature, detuning and potential structure. 

Without the partial waves, the magnesium situation 
would be an excellent application case for the two-state 
model. Now the question is, how well this approximation 
translates to strong fields and partial waves. 

2.5 Partial waves 

An example of the potential structure near the Condon 
point is shown in Fig. Ola). The couplings between the 
states form in fact a net of closely placed Condon points. 
The spatial separation between the individual two-state 
resonances suggests the possibility to partition the prob- 
lem into separate two-state processes, which is exactly the 
approach adopted in weak field studies ^H] ■ In this picture 
the approaching wave packet meets a sequence of two-state 
level crossings. It is then tempting to apply the Landau- 
Zener model to each crossing independently. 

In our specific magnesium case the separation into a 
sequence of two-state models would have a qualitatively 
important aspect. Because the Q branch is missing, there 
is a possibility for a " channelling effect" . Now each ground 
state partial wave (with the exception of the s-wave) is 
coupled to two excited rotational states, but at different 
Condon points. Assume that the coupling Vji is strong 
enough that Plz is almost unity for all crossings. Any part 
of the wave packet that arrives on a ground-state partial 
wave with I = I' > will thus move to the excited state 
j = I'—l at the first crossing. But in the next crossing, the 
excited rotational state is strongly coupled to the ground 
state partial wave with l = j—l = l' — 2 (see Fig.|2J). Thus 
the wave packet returns quickly to the ground state, where 
it climbs the centrifugal barrier, gets reflected, and follows 
the same path / — I' — 2 — ► j = I' — 1 — * I = I' backwards. 
This means that the collision is basically reduced to an 
elastic one for all partial waves except the s-wave. 



(a) 




ZkFUX 

Fig. 3. (a) Ground and excited state potentials for partial 
waves/rotational states up to I = 20 and j = 21. Open circles 
indicate the point of couplings forming a net of Condon points, 
(b) The corresponding dressed potentials. The parameters are 
S = -3-Tat and n = ir at . 



At this point we note that in our simulations no such 
strong suppression was seen to occur. Basically, the two- 
state approximation fails because even at moderate fields 
the other states must be taken into account. This is best 
seen if one considers the dressed states, i.e., potential 
curves that correspond to the eigenstates of the coupled 
quasimolecule+light system, see Fig. G^b). In order for 
the separation to independent two-state crossings to work 
out, the dressed states should show avoided crossings near 
each Rc- But it turns out that, for this to happen, the cou- 
pling strengths should be so low that Plz is far from unity. 
Therefore the division into separated, two-state Landau- 
Zener crossings is not applicable to the problem at hand, 
except for giving a rough tool to define the concept of 
strong field. 

Strong fields correspond to the adiabatic limit, which 
means that the wave packet should follow the dressed 
states smoothly, and position-dependent unitary transfor- 
mations describe the change from the bare state basis to 
the dressed state basis and back. The concept of adiabatic 
following and the basis change symmetry is broken by the 
spontaneous emission. It is still possible to treat the prob- 
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lem in the dressed state basis, but then the decay can go 
both ways as now practically all dressed states are super- 
positions of the ground and excited states. In the Monte 
Carlo method this leads to a complicated description, and 
therefore the bare states are usually preferred, even at 
strong fields. 

In the near-resonant situation the dressed and the field- 
free (bare) quasimolecule states do not match asymptot- 
ically at R — > oo, which makes the setting of the initial 
conditions difficult. If we set the wave packet initially on 
one of the dressed states, which kinetic energy should we 
give it, the initial ground state one, or should we take into 
account the energy difference between the bare ground 
state, and the selected dressed state? Spontaneous emis- 
sion allows us to circumvent this problem. 

Since the atoms approach slowly, the quasimolecule 
reaches a steady state (in terms of the ground and excited 
state populations) before the wave packet reaches the Con- 
don point. The approximate time scale for steady state 
formation is a few times fi/ T a ±. It does not really matter 
whether we set the system on the bare or dressed state ini- 
tially as it quickly adjusts itself to the steady state. Thus 
we can set the initial wave packet on a bare ground state 
corresponding to some particular value of I. It should be 
noted that in the Monte Carlo approach, the steady state 
appears only in the ensemble average. 

For practical reasons, we also need to truncate the 
infinite sets of partial waves and rotational states into 
finite ones, with some maximum values. Because of the 
centrifugal potential, only low-lying partial waves should 
contribute to a collision event. The maximum classically 
allowed angular momentum is obtained by matching the 
initial kinetic energy Ek to the centrifugal potential at the 
Condon point Rq (taken to be the s-wave Condon point, 
for simplicity), resulting in the value 



^max, class — 

Excitation can occur also for higher partial waves than 
this classical value, because of the quantum nature of the 
dynamics. The wave packet has a finite width, i.e., there 
are velocity components for which E > Ek- Also, there ex- 
ists a kind of tunnelling effect for atoms that nearly reach 
Rc, which is quickly suppressed as I increases. To be on 
the safe side, we have chosen l max to be a few partial waves 
higher than the classical value; this has been found to be 
a satisfactory description in previous work |2(j| , where the 
role of the partial wave structure was studied with wave 
packets assuming only independent Q branch (I — ► j = I) 
processes. 

It is important to avoid unwanted artificial effects that 
may arise because of the partial wave/rotational state 
truncation. In particular, choosing an odd number of par- 
tial waves may lead to the existence of a dark eigenstate 
in the Hamiltonian Such a state is an almost pure su- 
perposition of ground-state levels and any population that 
occupies it will survive, unaffected by decay, and is not ex- 
cited at the Condon point. If the system has a dark state, 
the frequent quantum jumps drive the quasimolecule pop- 
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Fig. 4. An example of the time evolution of the average ki- 
netic energy in a single collision simulation (an ensemble of 64 
realizations). The parameters are: S — — 3.0F a t, O = 1.0F a t, 
hi = 10.0fc r , Z max = 10 and the initial component is h — 8. 

ulation efficiently into it well before the actual collision 
occurs. 

Since a dark state can only exist if the Hamiltonian 
consists of an odd number of coupled levels, any suppres- 
sion of radiative heating because of a dark state will be a 
pure artifact of the finite matrix and must be avoided. In- 
deed, our numerical computations indicate that in the case 
of few partial waves, the excited-state flux is much sup- 
pressed if the number of states is odd, but the difference 
between the odd and even cases is less pronounced when 
the number of levels is larger. All the results reported here 
are obtained using an even number of states. 



3 Heating models 
3.1 Single-collision model 

In order to obtain information on the heating rate we cal- 
culate the average kinetic energy increase due to radiative 
heating in a single collision (see Appendix A.l for further 
details) . In this model the initial wave packet is reasonably 
well defined in position and momentum, so that we can as- 
sign the initial relative velocity v with the mean velocity of 
the packet, and the atoms are clearly separated. Next the 
atoms collide and heat up due to the radiative heating pro- 
cess, and move apart again. The single-collision simulation 
is stopped when the collision region is emptied of popula- 
tion and before the quasimolecule wave packet reaches the 
edge of the simulation space. From an ensemble of such 
simulations we obtain the average kinetic energy increase 
per collision AE SC . 

This approach is conditional to the point that the post- 
collision behavior of the wave packet follows the above 
classical picture; an important aspect of our simulations 
is to show that this is indeed the case. Figure 0] shows 
an example of the time evolution of the average kinetic 
energy and we see a clear stepwise behavior of kinetic en- 
ergy. Then AE SC can be calculated as the height of the 
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t(r') 

Fig. 5. An example of the time evolution of the mean position 
of the wave packet, (R), in a simulation box of L — 4A (solid 
line), and L = 2 A (dashed line). The parameters are: 8 — 
— 3.0-Tat, fi = 1.0-Tat, ki = lO.Ofc,-, 'max = 10 and the initial 
component is U = 8. 

energy step. The rough following of a classical trajectory 
is also demonstrated by the mean value of the wave packet 
position R (solid line in Fig. 0. Finally we note that the 
postcollisional distributions \&(R, t)\ 2 in all simulations 
clearly show that the collision region is emptied after the 
collision event. 

The inelastic process spreads the wave packet in mo- 
mentum and thus in position, and eventually the fast parts 
get reflected first at the edge of the simulation space (at 
R = 4A here), and (R) loses its classical-like behavior 
(Fig. |5| . Similarly, the wave packet is delocalized over the 
whole simulation space, and remains like that. This is re- 
flected in the long-time behavior of (R) (dashed line in 
Fig.[5J|. Clearly the classical- like trajectory description ap- 
plies only to the first collision. For subsequent collisions we 
can not define any classical-like asymptotical initial and 
final states. 

Results from single-collision calculations are related to 
the physical heating rate coefficient of a thermal cloud at 
temperature T through the relation (cf. Appendix A.l) 

£ (2k + l)AE ac (E k ,k), (11) 

l(even) 

where Qt = (Sir/iksT '//i 2 ) 3//2 is the translational partition 
function, and AE sc (E k , l{) is the numerical kinetic-energy 
increase calculated for the initial kinetic energy E k and 
initial ground state partial wave In practice we calculate 
the nonaveraged, energy-dependent rate coefficient 

K H (E k ) = -^- V {2l i + l)AE sc (E k ,l i ). (12) 

t(even) 

The multi-level Monte Carlo computations are extremely 
time consuming, and we must therefore be content to cal- 
culate the nonaveraged rate coefficients Ku(E k ) for only a 



few instances of parameter values. We argue that Kii{Ek) 
still indicates the order of magnitude of the heating rate 
for temperatures around T = Ek/ks- 

For the weakly-coupled case that has been treated pre- 
viously the wave packet populates the ground-state 
partial wave on which it was initially concentrated, h, and 
in the collision region it populates also the neighboring 
excited-state levels, but is not allowed to spread beyond 
those. However, as the Rabi coupling Q is increased, one 
quickly enters the regime where the population is dis- 
tributed among several ground-state partial waves before 
the collision. The dependence on the initial state is there- 
fore rather weak. For example, when Q = 1.0-Tat, the dif- 
ference in kinetic-energy increase between the cases l\ = 
and l{ — 8 is only 25% (see Table |2J ■ These dependencies 
are also affected by the choice of initial position, but the 
dependence of the final partial- wave-summed rates on the 
initial position appears to be weak. 

3.2 Multicollision model 

An alternative method to calculate the heating rate co- 
efficient is a multicollision simulation. The essence of the 
method is that the quasimolecule population is reflected 
at the edge of the numerical space and thus performs re- 
peated collisions. The beginning of the simulation resem- 
bles the single-collision situation (see Fig. EJ). But then 
multiple collisions begin to occur and the population is 
quite rapidly distributed over the whole simulation space. 
In a sense, the wave packet components slosh back and 
forth, spreading in momentum as they reach the Condon 
point, and spatial coherence is quickly lost. 

The quasimolecule gains kinetic energy when collisions 
occur, the wave packet moves faster, and the collision rate 
increases. We can now monitor the time evolution of the 
average kinetic energy and calculate the heating rate di- 
rectly from the curve. The important point here is that 
the memory effects from initial conditions are lost and the 
result of the simulation should be independent of them. 
The quantum jump processes contribute strongly in mix- 
ing the population distribution of various partial waves 
between collisions. This brings the multicollision model 
closer to the actual situation where atoms emerging from 
a collision on some particular partial wave meet individu- 
ally other atoms, forming a collision complex with partial 
waves that have nothing to do with the partial waves of 
the previous collision. 

In the very beginning of a multicollision simulation we 
have the single-collision stepwise behavior in the average 
kinetic energy, but after the spatial coherence (i.e., local- 
ization to the vicinity of a classical trajectory) is lost, we 
find that the kinetic energy increases linearly, see Fig. [fj] 
This allows one to obtain the heating rate as the slope of 
the curve for the particular density corresponding to the 
size of the simulation space, as described in Appendix A. 2. 
The increase in energy becomes non-linear later on. This is 
because the energy increase implies a faster collision rate 
and better survival, which compete with the weaker exci- 
tation probability (as discussed in Sec. 12.41) . But we have 
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found that at that point the atoms are already quite hot 
and strongly affected by Doppler cooling in any case. In 
order to obtain the heating rate in practice, we measure 
the slope of the energy in a time interval [ti,t2], where 
the initial time t\ is chosen so that the wave packet has 
already had time to rid itself of any dependence of ini- 
tial state, and the upper limit ti is chosen well before the 
curve becomes non-linear. 

The numerical result from the multicollision computa- 
tion is related to the actual heating rate coefficient simi- 
larly to Eq. lO, c f- Appendix A. 2: 



K H (E k ) = 



E k 
hQE, 



l(even) 



i)- 



dE mu i 2L 



dt 



(13) 



Here, L is the length of the simulation grid, v is the ini- 
tial velocity v = \j2Ekj [i, and dE mn \ is the energy of 
the wave packet in the multicollision model. Thus, the 
quantity under the summation sign is the energy gain per 
collision (the factor 2 multiplying L reflects that the wave 
has to travel across the grid twice between two collisions) . 
Moreover, we have already noted that in the multicollision 
model the results are more or less independent of the ini- 
tial partial wave; thus, we directly perform the sum over 
even L and obtain 



K H {E k ) = \{h 



c+1) 



, s/2% 2 2L dE, 



mul 



dt 



V22LRj dEmul 



dt 
(14) 

As an internal check, the two simulation methods should 
be compared by matching the energy gain per collision: 



AE SC = 



dE mu i 2L 
dt v 



(15) 



The left- and right-hand sides should be summed over par- 
tial waves before they can be compared, but anticipating 
the result that the single-collision simulations are almost 
independent of the initial angular momentum, the depen- 
dence can be neglected. 



4 Numerical results for heating rate 
coefficients 

The parameters for the simulations are chosen to corre- 
spond to realistic experiments on magnesium. The atomic 
linewidth is -T a t = 391-Er, and the recoil energy En = 
ks x 9.8 /iK. We have chosen to work with a detun- 
ing 8 between 1 and 3 times the atomic linewidth _r at , 
and varied the Rabi coupling ft between 0.1 and 2 times 
r a t, corresponding to laser intensities of 36 mW/cm 2 to 
14.2 W/cm 2 . The initial momentum ki has been chosen 
to be ki = lQk r in most instances, giving an energy E% = 
10QE R = k B x 0.98 mK. This is slightly below the Doppler 
temperature To = 1.9 mK. We present the results of each 
of the two simulation methods, and compare them both 
to each other, and to the outcome of a simpler two-state 
model. 
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Fig. 6. An example of the time evolution of average kinetic 
energy in the multicollision simulation. The kinetic energy be- 
haves linearly and the slope of the curve gives the heating rate. 
Here 5 = -3.0r at , O = 1.0r at , h = lO.OAv, ^max = 10 and the 
initial component is k = 0. 



4.1 Average kinetic energy increase 

The results for the average kinetic energy increase in the 
single-collision model are displayed in Table [3 At the top 
of the table are the results from using only two partial 
waves, the ground-state partial wave as indicated in the 
column l{ and the neighboring higher excited-state level, 
j = k + 1. Below the horizontal bar are the results from 
the full calculation, incorporating the first 6 or 7 partial 
waves (i.e., up to i max = 10 or 12). The initial mean po- 
sition is set to Ro = 5\/2ir. Some clear differences can 
be seen. First, when starting on the s-wave, the energy 
increase is larger in the two-state case than in the mul- 
tistate case. The reason is rather trivial: in the two-state 
case there is only one channel for excitation; the probabil- 
ity for this transition is proportional to 2/3, as can be seen 
in Table rj] On the other hand, in the multistate calcula- 
tion, many partial waves are occupied as the wave packet 
enters the crossing region, and the transitions to the ex- 
cited state take place with a suppressed probability as can 
also be seen in Tabled when I increases towards infinity, 
a 2 ; approaches 1/6. Conversely, when k is high enough, 
the multistate result for the energy increase is higher be- 
cause the population now spreads to lower partial waves 
before the collision. Figure0illustrates this. These results 
demonstrate that in the multistate case the redistribution 
of population before the collision leads to a weak depen- 
dence on the initial partial wave, as already discussed in 

sec.nm 

The intensity dependence of the multistate single-colli- 
sion results is demonstrated in Table |31 and we discuss it 
later, in Sec. 14.31 



4.2 Multicollision heating rate 

For the multicollision simulations, the size of the simula- 
tion box is set to L — 2. OX. Again, a two-state calculation 
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Fig. 7. Energy increase due to radiative heating in the single- 
collision model as a function of initial partial wave. Crosses 
indicate the two-state result and circles are the results from 
using the full partial- wave structure (multistate model). The 
data refers to Tabled 

Table 2. The comparison of results from two-state and mul- 
tistate single-collision models. Here Ni is the total number of 
states in the simulation (either 2 or 12), k is the initial state, 
ki the initial momentum, and AE sc (Efi) the average increase 
of kinetic energy. The energy increase is calculated as a time 
average of the kinetic energy in the region where it is flat, sub- 
tracted by the initial energy. The error refers to statistical error 
of the kinetic energy at the end of the Monte Carlo simulation. 





fi(r at ) 


<5(r at ) 


k 


ki {kr ) 


AE SC (E R ) 


2 


1.0 


-3.0 





-10 


705 ± 106 


2 


1.0 


-3.0 


6 


-10 


337 ± 44 


2 


1.0 


-3.0 


8 


-10 


371 ± 40 


2 


1.0 


-3.0 


12 


-10 


187 ± 33 


12 


1.0 


-3.0 





-10 


404 ± 51 


12 


1.0 


-3.0 


8 


-10 


307 ± 31 


14 


1.0 


-3.0 


12 


-10 


332 ± 33 



Table 3. Intensity dependence of the results from the mul- 
tistate single-collision simulations. The parameters are as in 
Table H 



N 


fi(r at ) 


d(r«o 


li ki {kr ) 


AE SC (E R ) 


12 


1.0 


-3.0 


8 -10 


307 ± 31 


12 


0.5 


-3.0 


8 -10 


145 ± 16 


12 


0.1 


-3.0 


8 -10 


16 ± 3 



has been done in order to compare with the full partial- 
wave structure up to Z max = 10. The results for the heating 
rate are displayed in Table 

Comparing the two- and twelve state results, there is 
a similar difference, as already seen in the single-collision 
results. Moreover, and like discussed earlier, the single col- 
lision results indicated only the weak dependance of the 
initial state. We expect this to hold also especially here for 
the multicollision model. Assuming that the initial state 
dependence is weak to effect the first collision, the role of 
the initial state should be even weaker for the following 



Table 4. Comparison of the two-state and multistate results 
in the multicollision model. Here h is the initial state, ki the ini- 
tial momentum, and dEm U i/dt the slope of the energy- versus- 
time curve. 





o(r at ) 




li 


ki {k r ~) 


d£ mul /dt(Sflr at //i) 


2 


1.0 


-3.0 





-10 


4.48 


2 


0.5 


-3.0 





-10 


2.23 


2 


0.4 


-3.0 





-10 


1.00 


2 


0.3 


-3.0 





-10 


0.50 


2 


0.2 


-3.0 





-10 


0.22 


2 


0.1 


-3.0 





-10 


0.069 


12 


2.0 


-3.0 





-10 


1.20 


12 


1.0 


-3.0 





-10 


0.79 


12 


1.0 


-3.0 


6 


-10 


0.67 


12 


0.5 


-3.0 





-10 


0.40 


12 


0.1 


-3.0 





-10 


0.039 



collisions which actually determine the slope of the heat- 
ing curve, i.e. the heating rate. This view is confirmed by 
the simulation result with li = 6 in Table 0| Thus the re- 
sults concentrate here on the intensity dependence of the 
heating rate which is discussed in the following subsection. 



4.3 Comparison of results 

The outcome of the single-collision and multicollision mod- 
els are now to be mapped to physical quantities using 
Eqs. lfH?|) and ijTlj l. Equation ijl"%|) is rewritten as 



K H {E k ) 



2ttE 



1/2 ' 



1 



.3/2 



AE sc (R.n.) 
fci(R.u.) : 



(16) 



where (R.u.) means recoil units. The fraction within brack- 
ets has the numerical value 1.1307x 10~ 43 Wm 3 . The heat- 
ing rates from the multicollision computations are con- 
verted into energy increase per collision with the conver- 
sion factor 2L/v = 491 in recoil units. The results are 
summarized in Table [3] The Monte Carlo computations 
are extremely time consuming, but because of the weak 
dependence on initial partial wave and the large error 
bars, it has not been necessary to perform the calcula- 
tion for each l[, but the values tabulated in Tables [21 and 0] 
have been sufficient. The correspondence between the two 
models is seen to be fair. 

If we assume fe; = 10k r , S = —3.01^, Vc = 2i?/3, take 
the leading i?~ 3 -term of the potential U(R) and ignore the 
partial waves (I = j = 0), the argument in the exponen- 
tial in the Landau-Zener term in Eq. @ becomes equal 
to 7.2 x (Q/r &t ) 2 . For Vc = ^/6, the asymptotic limit for 
the aji for large I, this decreases to 1.8 x (Q/r at ) 2 . This 
means roughly that in simulations with J? = _T at we are 
more or less saturating the Landau-Zener two-state result 
for all partial waves. Although the model is not directly 
applicable to the multistate case, this result indicates nev- 
ertheless, that we should be in a region where the usual 
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Table 5. Physical heating rates in units of 10 _40 Wm 3 given by 
the two numerical methods. R s stands for single collision and 
Rmui for multicollision heating rates. All the results are com- 
puted in the twelve-state model allowing for the full partial- 
wave structure. 



f2(r at ) 


5(r at ) 


E k /k B (mK) 


R s 




2.0 


-3.0 


0.98 


3.74 


4.03 


1.0 


-3.0 


0.98 


2.10 


2.65 


0.5 


-3.0 


0.98 


0.99 


1.34 


0.1 


-3.0 


0.98 


0.11 


0.13 


1.0 


-3.0 


3.92 


1.21 




0.5 


-1.0 


0.98 




0.23 



weak-field expectation of linear increase of the excitation 
rate with intensity / fails. Figure |S] plots the heating rate 
as a function of the Rabi coupling for 5 = — 3-T a t and 
ki = — 10k r . We see an increase which is at first roughly 
linear in fl, but then becomes even weaker, indicating a 
possible saturation. This is in qualitative agreement with 
two-state models [TJ; when the Landau-Zener process sat- 
urates, the re-excitation of decayed population leads to an 
increase in excited state survival, which is linear in fl j.'iOj . 

5 Discussion and conclusions 

We have explored in detail the regime of strong laser fields 
in cold alkaline-earth collisions, which necessitates taking 
the full partial- wave structure into account. For Rabi cou- 
plings fl of the order of or larger than the linewidth _T at , 
the wave function of the colliding atoms quickly spreads 
and populates all partial waves before and during the col- 
lision, yielding results that depend only weakly on initial 
conditions. 

The outcome of the two computational models, the 
ordinary simulation of a single collision and the multicol- 
lision model where the wave sloshes back and forth in a 
finite space, are seen to match well. A priori one tends to 
associate the multicollision model as a more rigorous way 
to calculate the heating rate. This is because consecutive 
collisions, during the system's long term time evolution, 
keep track on the changes of the system properties in a 
more realistic way compared to the single collision model. 
However, since the results produced by the two methods 
have reasonable agreement, we come to the conclusion 
that the computationally lighter single collision method 
provides a more practical approach. 

It is worth noting that recent experiments on alkaline- 
earth atoms in magneto-optical traps show that the Dop- 
pler limit is not reached 3lE3SSj- Our earlier results indi- 
cate that the radiative heating does not prevent reaching 
of the Doppler limit for current rather low atomic densities 
in MOTs The improved method for cold collisions, 
multicollision model presented here, gives further confir- 
mation for this view. The authors of Ref. attribute the 
extra heat to the transverse spatial intensity fluctuations 
of their one dimensional optical molasses. However, a full 
three dimensional set-up may introduce additional factors 
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Fig. 8. Numerically calculated energy-dependent heating rate 
R as a function of Rabi coupling O for the energy E k — ks x 
0.98 mK and the detuning 5 — — 3-T at . Crosses are the results 
of the single-collision simulations and circles are the results of 
the multicollision model. 



to the problem. From the point of view of atomic collisions, 
the partial wave treatment includes three dimensional as- 
pects, and the conclusions of this paper remain the same. 

In conclusion, we have developed methods, based on 
Monte Carlo wave function simulations, to include partial 
waves to the cold collision problems in strong laser fields. 
The methods have been described in detail highlighting 
the most relevant aspects of the problem. For example, 
the appropriate number of partial waves has to be even 
to avoid the artificial dark state to influence the dynam- 
ics. The results show the linear increase of the heating 
rate with laser intensity, even in the regime where the 
Landau-Zener excitation probability saturates. This is in 
qualitative agreement with two-state models when the re- 
excitation of decayed population has been accounted for. 
We think that the work presented here presents a con- 
siderable step towards a complete partial wave treatment 
of cold collisions in the presence of near-resonant light. 
This topic is also relevant in the general framework of fun- 
damental collision theory. A completely satisfactory and 
realistic model remains as a challenging task for future 
studies. 
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Appendix A: Rate coefficients 
A.l Single collision rate 

The collisional heating rate Kjj(n) describes the kinetic 
energy change per unit volume and unit time. For inden- 
tical particles it will be equal to Knn 2 /2, where Kh is 
the rate coefficient for heating, and n is the atomic den- 
sity. The factor of 1/2 removes the doubling in collision 
counting for identical particles. The total collisional heat- 
ing rate is obtained by integrating kh over trap volume. 
For simplicity we ignore the trap loss, so n is constant in 
time. In any case, since the rate coefficient Kh is density- 
independent, any time dependence in n will not affect it. 

For trap loss rate coefficient K\ oss one normally cal- 
culates the collision frequency times the probability for 
loss, i.e., dctioss, where v is the relative velocity, and <7i oss 
is the cross-section for collisions, including the probability 
for a loss event to occur. This simple classical picture is 
connected to thermodynamics by assuming a distribution 
of velocities, f(v)dv, over which we take an average, and 
obtain K\ oss = (fcioss)' The dynamics of the two-body col- 
lision on the microscopic level then enters in calculating 

^loss ■ 

To calculate the heating rate Kh is slightly more com- 
plicated, because in addition to the probability of an in- 
elastic collision event to happen, we also need to estimate 
the amount of kinetic energy increase associated with it. 
Technically, we would have a continuous distribution of 
final energies corresponding to each value of v, due to the 
randomness of the spontaneous emission events. For prac- 
tical reasons we consider an averaged rate, i.e., to calculate 
the average change in kinetic energy per unit time. Since 
we can perform the averaging over final energy states be- 
fore the averaging over initial states, we define a heating 
cross section an(v) (units energy x distance 2 ), that gives 
the difference of the average final relative kinetic energy 
and initial relative kinetic energy. 

In the partial wave approximation we can write the 
quantum mechanical cross-section for identical atoms in a 
three-dimensional trap as 12 



oo 




(17) 



where k is the wave number related to v and Piiv) is the 
event probability. Thus, as a generalization, we write for 
heating 



where ks is the Boltzmann coefficient, /i is the reduced 
two-particle mass (equal to m a t/2 for identical atoms), 
and T is the gas temperature at equilibrium. By defining 
a dimensionless energy ratio e = \jJLV 2 / '(ksT) we get 



7 rp />OG °° 

kQt Jo n 



(20) 



If we truncate the partial wave sum e.g. with the clas- 
sical trajectory argument, then the upper limit for the 
partial wave sum becomes energy dependent. Therefore it 
is preferable to keep the sum inside the energy integral. 

Due to the complexity of the Monte Carlo simulations 
we are, in practice, limited to calculating AE sc (v, I) for a 
rather narrow range of initial relative velocities v. As in 
Ref. rjHj ! w e define a single-energy rate as 



K H ( S ) = ^J2(2l + l)AE sc (e,l). 



(21) 



As one can see immediately, this mimics the case where 
e = l, and the energy dependence on the integrand comes 
only from the Boltzmann factor. 



A. 2 Multicollision rate 

The classical picture in the derivation of a rate coefficient 
is a particle with velocity v, travelling a distance vAt in 
time interval of length At, and having the "area" given by 
the cross section a. With gas density n, it meets v Atari 
particles in that time, so that the collision frequency for 
one particle becomes van. The density of such particles 
brings in another n. 

In the multicollision model the fixed simulation box 
size L gives an alternative way to approach the problem. 
Now the effective collisional volume is 2La, in which the 
multiple collisions eventually become inseparable, and the 
kinetic energy does not increase in a stepwise manner. 
Here the quantum jumps and steady state formation make 
sure that apart from the change in the kinetic energy, there 
is no memory effects between collisions. Therefore, on av- 
erage the situation corresponds to multiple inelastic colli- 
sions when we have a particle density l/(2La), although 
in reality the particles do not collide at regular intervals. 
If we normalise the rate of change to density, we get 



K H {e) 



k B T f 



hQ: 



dE2L 
dt v 



(22) 



^HpfE^+^-CM), (is) 

o 

where AE sc (v,l) is the average single-collision energy in- 
crease related to the "initial" partial wave 

In an isotropic Maxwcll-Boltzmann gas the distribu- 
tion of relative velocities v is given by 

f(v)dv = 4tt (-JL-^ /2 Je-^l^dv, (19) 



There is one caveat, though. If the change in kinetic 
energy is a strong function of v, then clearly for suitably 
long times ^ is not a constant, and one would need to 
revert to calculating heating rates from single-collision re- 
sults. The fact that in simulations ^M- remains a constant 

at 

for a wide energy region makes it possible to use results 
from multicollision studies in the manner descibed above. 
More importantly, it implies that as one approaches the 
Doppler cooling temperature, the heating rate per partial 
wave remains constant. 



12 



J. Piilo et al.: Cold collisions in strong laser fields: partial wave analysis of magnesium collisions 



References 

1. H. J. Metcalf and P. van der Straten, J. Opt. Soc. Am. B 
20, 887 (2003). 

2. M. H. Anderson, J. H. Ensher, M. R. Mattews, C. E. Wie- 
man, and E. A. Cornell, Science 269, 198 (1995); K. B. 
Davis, M.-O. Mewes, M. R. Andrews, N. J. van Druten, 
D. S. Durfee, D. M. Kurn, and W. Ketterle, Phys. Rev. 
Lett. 75, 3969 (1995); C. C. Bradley, C. A. Sackett, J. J. 
Tollett, and R. G. Hulet, Phys. Rev. Lett. 75, 1687 (1995). 

3. B. DeMarco and D. S. Jin, Science 285, 1703 (1999). 

4. A. G. Truscott, K. E. Strecker, W. I. McAlexander, G. 
Partridge, and R. G. Hulet, Science 291, 2570 (2001). 

5. K.-A. Suominen, J. Phys. B 29, 5981 (1996). 

6. J. Weiner, V. S. Bagnato, S. Zilio, and P. S. Julienne, Rev. 
Mod. Phys. 71, 1 (1999). 

7. C. J. Pethick and H. Smith, Bose-Einstein Condensation 
in Dilute Gases (Cambridge University Press, Cambridge, 
2001). 

8. The scattering rate of photons sets the Doppler limit, 
which is 2 mK for the strong cooling transition So- Pi 
in magnesium, in comparison with the photon recoil limit, 
10 fiK. In alkali atoms the degeneracy of the ground state 
allows one to break the Doppler limit and reach the recoil 
limit with tools like Sisyphus cooling, but this will not work 
for the alkaline earth atoms that have no such degeneracy. 

9. F. Bardou, J. -P. Bouchard, A. Aspect, and C. Cohen- 
Tannoudji, Levy Statistics and Laser Cooling (Cambridge 
Univ. Press, Cambridge, 2002). 

10. M. Kasevich and S. Chu, Phys. Rev. Lett. 69, 1741 (1992). 

11. A. Gallagher and D. E. Pritchard, Phys. Rev. Lett. 63, 
957 (1989). 

12. P. S. Julienne and F. Mies, J. Opt. Soc. Am. B 6, 2257 
(1989). 

13. P. S. Julienne and J. Vigue, Phys. Rev. A 44, 4464 (1991). 

14. K.-A. Suominen, Y. B. Band, I. Tuvi, K. Burnett, and P. 
S. Julienne, Phys. Rev. A 57, 3724 (1998). 

15. G. Hillenbrand, C. J. Foot, and K. Burnett, Phys. Rev. A 
50, 1479 (1994). 

16. B. M. Garraway and K.-A. Suominen, Rep. Prog. Phys. 
58, 365 (1995); Contemp. Phys 43, 97 (2002). 

17. H. J. Carmichael, Statistical Methods in Quantum Optics 
Vol. I, (Springer, Berlin, 2002). 

18. J. Dalibard, Y. Castin, and K. M0lmer, Phys. Rev. Lett. 
68, 580 (1992); K. M0lmer, Y. Castin, and J. Dalibard, J. 
Opt. Soc. Am. B 10, 524 (1993). 

19. M. Machholm, P. S. Julienne, and K.-A. Suominen, Phys. 
Rev. A 64, 033425 (2001). 

20. K. Sengstock, U. Sterr, J. H. Miiller, V. Rieger, D. Better- 
mann, and W. Ertmer, Appl. Phys. B 59, 99 (1994). 

21. D. N. Madsen and J. W. Thomsen, J. Phys. B. 35, 2173 
(2002). 

22. J. Piilo, E. Lundh, and K.-A. Suominen, Phys. Rev. A 70, 
013410 (2004). 

23. W. J. Meath, J. Chem. Phys. 48, 227 (1968). 

24. W. J. Stevens and M. Krauss, J. Chem. Phys. 67, 1977 
(1977). 

25. E. Czuchaj, M. Krosnicki, and H. Stoll, Theor. Chem. Acc. 
107, 27 (2001). 

26. M. Holland, K.-A. Suominen, and K. Burnett, Phys. Rev. 
Lett. 72, 2367 (1994); Phys. Rev. A 50, 1513 (1994). 



27. In fact in simulations we also determine from which partic- 
ular state j the collapse originated, and along which branch 
it took place, so after a jump only one ground state partial 
wave is occupied. 

28. J. Piilo, K.-A. Suominen, and K. Berg-S0rensen, J. Phys. 
B 34, L231 (2001); Phys. Rev. A 65, 033411 (2002). 

29. J. Piilo and K.-A. Suominen, Phys. Rev. A 66, 013401 
(2002). 

30. When the Landau-Zener process saturates and re- 
excitation dominates, we can assume the delayed de- 
cay model |14|. where the exponential decay starts only 
once the quasimolecule has reached the edge of the re- 
excitation zone. The location of this edge is roughly where 
\U(R) — U(Rc)\ — fi, so the energy increase obtained by 
the quasimolecule before the exponential decay is propor- 
tional to Q. 

31. X. Xu, T. H. Loftus, M. J. Smith, J. L. Hall, A. Gallagher, 
and J. Ye, Phys. Rev. A 66, 01 1401 (R) (2002). 

32. FY. Loo, A. Brusch, S. Sauge, M. Allegrini, E. Arimondo, 
N. Andersen, and J. W. Thomsen, J. Opt. B: Quantum 
Semiclass. Opt. 6, 81 (2004). 

33. T. Chaneliere, J.-L. Meunier, R. Kaiser, C. Miniatura, and 
D. Wilkowski, J. Opt. Soc. Am. B 22, 1819 (2005). 



